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We investigate the scaling behavior of the cluster size distribution in the Drossel-Schwabl Forest 
Fire model (DS-FFM) by means of large scale numerical simulations, partly on (massively) parallel 
machines. It turns out that simple scaling is clearly violated, as already pointed out by Grassberger 
[P. Grassberger, J. Phys. A: Math. Gen. 26, 2081 (1993)], but largely ignored in the literature. 
Most surprisingly the statistics not seems to be described by a universal scaling function, and the 
scale of the physically relevant region seems to be a constant. Our results strongly suggest that the 
DS-FFM is not critical in the sense of being free of characteristic scales. 



1. INTRODUCTION 

The Drossel-Schwabl Forest Fire Model ^ (DS-FFM) 
is one of the paradigms of non-conservative SOC Its 
importance comes primarily from the fact that the model 
has non-conservative microdynamics. It therefore an- 
swered the question whether conservation is necessary 
for criticality in driven systems. 

The claim that the DS-FFM is critical comes from the 
fact that it shows powerlaw-like behavior for several ge- 
ometrical properties of the dissipation events. The aver- 
age size of these is divergent in the so-called SOC limit, 
where all timescales get separated so that the rate of the 
external drive becomes infinitely slow and the total in- 
flow diverges. If one assumes stationarity this is trivial 
to prove jl| H, 1^. However, as usual in numerical simu- 
lations, it is not possible to investigate the model in the 
limit of divergent drive (0 ^ in the notation below), 
as finite size limits the correlation length and therefore 
destroys any possible criticality . It is remarkable that 
most of the literature available for this model is mainly 
concerned with finding critical exponents and identify- 
ing supposedly critical quantities. It seems that no au- 
thors question whether the model is critical at all and 
if so in which sense. In this paper we carefully investi- 
gate the "scaling function" of the cluster size distribution 
and show that it is indeed an open question whether the 
model is truly critical: Not only is there no way to prove 
its criticality, there is also numerical evidence that the 
model does not become scale free. 



2. DEFINITION OF THE MODEL AND 
METHODS 

The model has been described several times and in 
great detail elsewhere [|l|, |3[ Q . Therefore, the description 
presented here is rather succinct. The model is defined 
on a d-dimensional lattice of linear size L, where each 
lattice site can be in one of two states, "occupied" or 
"empty". The lattice is then updated according to the 
following rules: 

• Driving: Choose randomly (1/6*) sites, one after 



the other. If its state is "empty" turn it into "oc- 
cupied" . 

• Relaxation: Choose one site at random. If it is 
empty continue with the first step. Otherwise make 
"empty" all the sites in the cluster the site chosen 
belongs to. In this case the update is considered to 
be successful. Continue with the first step. 

Here a cluster is defined in the usual fashion as the set 
of occupied sites which are connected via nearest neigh- 
bor interactions, i.e. two sites belong to the same cluster 
if they are nearest neighbors or there is a path between 
them along sites, which belong to the same cluster. We 
have applied periodic boundaries in all our simulations 
and restrict ourselves to the two dimensional square lat- 
tice. 

The cluster removed in the second step is called the 
"burnt cluster". To measure the overall distribution of 
clusters within the system, one usually measures the size 
of the burnt cluster Q , the distribution of which is bi- 
ased by a factor s. To see that, we define n(s) to be the 
ensemble-averaged, site-normalized density of clusters of 
size s in the system. Then, the probability that a ran- 
domly chosen site is connected to a cluster of size s is 
sn{s), as in standard percolation This distribution is 
probably the most important in the model. Other quan- 
tities are the distribution of the burning time, which is 
defined as the maximum Manhattan distance (shortest 
path on the square lattice) from the initially chosen site 
of a burnt cluster to all other sites in the same cluster, 
and the correlation functions as defined and discussed in 
. In this paper we solely concentrate on the distribu- 
tion n{s). 

Using a new implementation of the model 0| we are 
able to simulate the system on very large scales and at the 
same time keep track of the entire distribution n(s), in- 
stead of measuring the biased distribution sn{s), as done 
usually [||. Between two updates the changes in n{s) are 
only of the order {1/0), so it is a highly correlated quan- 
tity. However, by using standard cluster labeling tech- 
niques |8|, it is possible to calculate the full histogram 
n{s) essentially without increasing the computing time, 
which depends almost exclusively on 9 and is essentially 
independent of the system size. Compared to the stan- 
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dard simulation, we gain up to two orders of magnitude 
in performance A similar method was recently in- 
troduced for standard percolation Using the same 
amount of computing time the results are significantly 
less noisy than those of the standard implementation (for 
example jl^, which we have used as reference to check 
the validity of our results). Large system sizes enable 
us to rule out any finite size effects as described below. 
The results have been partly cross-checked using a differ- 
ent random number generator (all results presented here 
make use of ran2 from pl] |, and for checking rani from 
has been used). 
Finite size effects have been ruled out by the follow- 
ing direct method: For each value of {I/O) a significantly 
larger system was simulated with exactly the same value 
for (1/0). The linear size L was typically increased by a 
factor 2. The smallest systems we used were L — 1000, 
the largest L = 32000. By comparing the histograms of 
different system sizes in conjunction with the standard 
deviation calculated for them, it is possible to decide 
whether a system size is affected by finite size effects or 
not. Compared to other simulations published p^ , 
which also claim not to suffer from finite size effects, our 
system sizes are huge. This simple method of "redoing" 
all simulations and using lattices which are much larger 
than actually needed has the obvious disadvantage of be- 
ing inefficient, but there is probably no more a direct way 
of identifying finite size effects ||lj] . This waste of com- 
puting power is overcompensated by the efficiency of the 
algorithm and self-averaging [fsl. 



3. RESULTS 

The focus of this paper is the presumably universal 
scaling function of the distribution n(s). Similar to the 
correlation function one expects 



(1) 



if "simple scaling" applies, which is already known not to 
be the case in the presence of finite size effects . The 
L-dependence of this quantity is omitted in the following 
wherever the context allows it. It is worthwhile to note 
that this is usually the definition of the exponent r. The 
function Q is the (presumably) universal scaling function, 
which depends only on the ratio s/sq{9), where sq(9) is 
the only scale of the distribution. This scale depends only 
on the external parameters, in our case 9. Simple scaling 
allows another scale, namely the lower cutoff, but this is 
fixed or at least bounded. 

The function g{x) is usually smooth for small values 
of X, therefore it does not make a big difference whether 
we investigate n(s) or 



n{s;9) = 



n{s;9) 
g{s/so{9)) 



(2) 



That this particular choice of the normalization does not 
affect the overall results can be seen in Tab. |, where 
the absolute value of n(l, 9) is listed for different values 
of 9. Also shown in this table is the first moment of 
the distribution or the average density of trees, which is 
defined as 



P - 



sn{s; 9) 



and the second moment of the distribution. 



J2s'^n{s;9) 



(3) 



(4) 



which is the average size of the cluster connected to a 
randomly chosen occupied site. 

Before presenting the actual results, we first discuss 
the numerical quality of the results. 



3.1. Avoiding finite size effects 
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FIG. 1: Ratio r{s; 9, Li, L2) = n{s;0, Li)/n{s- 9, L2) with 
(1/9) — 2000 for two pairs L\,L2 with error-bars (one stan- 
dard deviation; the error-bars as well as the data shown are 
exponentially binned). The data are from short runs (10^ up- 
dates for statistics). Finite size effects have been considered 
negligible under the condition that (almost all) error-bars for 
this ratio have covered 1 (marked by a dashed line) in the 
relevant range, a) Li = 4000 and L2 ~ 8000: Almost no 
finite size effects, the deviation from 1 is probably due to 
noise. Note the fine scale of the ordinate, b) Li = 1000 and 
1/2 = 8000: Systematic, strong finite size effects for s ^ 10'*. 
The scale of the ordinate is five times larger than in a). Data 
of this quality have been dismissed. 

Throughout this paper we initially performed 5 • 10^ 
successful updates (as defined in Sec. 2) as transient (and 
therefore rejected them) and the same number for pro- 
ducing statistics, apart from runs for calculating error- 
bars, where only 10^ updates has been used for statis- 
tics, see below. It is known that the transient can be 
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TABLE I: Static quantities for different choices of L and 
(1/S). The estimation of the standard deviation of the tree 
density p, cr^(p) = (p^) — (p)^, where the average runs over 
the ensemble, is unfortunately based only on a small subset 
of the configurations produced, and in the case of the large 
systems {L > 16000) only on a fraction of the lattice. How- 
ever, it is apparent that it behaves like 1/L, as expected for 
a system without finite size effects. The density of clusters 
of size 1, n(l), serves as the normalization of n. The average 
cluster size is denoted by (s), for definition see but due 
to a truncation in the histogram for some of the simulations 
in the range 2000 < (1/9) < 16000, the number presented 
is actually the average size of the burnt cluster. In the sta- 
tionary state it is - apart from statistical fluctuations - also 
given by (1 - {p))/{e{p)) §. Values of (1/6) and L printed 
in bold indicate results shown in Fig. ^, the other results are 
only for comparison. All data are based on 5 • 10^ (successful) 
updates (s. Sec. 2) for the transient and statistics, apart from 
those printed in italics which are based on short runs (5 • lO'' 
updates for the transient and 1 • 10*' updates for statistics). 
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FIG. 2: The binned histogram n{s;0,L) for two different 
values of L and fixed 6 as in Fig. |l^. In this plot the two his- 
tograms are virtually indistinguishable. However, note that 
the deviations shown on Fig. nib would also hardly be visible 
in this type of plot, as shown in the inset. 



very long (note that the time unit in Q is expressed 
in our units by multiplying it with (\ / 9) / {pL"^)) , but in 
all cases presented the number of initial steps seemed to 
be more than sufficient. Numerical checks indicate that 
the cluster size distribution is very stable against the size 
of the transient, i.e. even a transient, which is presum- 
ably too short, still produces reasonable results for n(s). 

All systems have been initialized by a random inde- 
pendent distribution of trees with density 0.41. 

The standard deviation of the binned histogram is not 
completely trivial to calculate. In particular, its compu- 
tation requires a significant amount of CPU time, and 
was therefore only calculated for the smaller system sizes 
(up to L = 8000) and in shorter runs (only 10^ updates 
for statistics, but 5 • 10^ for transient). We resorted to 
visual examination for the larger systems when compar- 
ing n{s;9,L) for different system sizes. Fig. pi and b 
show the ratio of n(s; 9, L) for two different system sizes. 
A deviation of this ratio from 1 indicates a difference 
in the statistics and therefore the presence of finite size 
effects. Fig. |l|a shows a typical case we accepted as rea- 
sonable agreement. Here Li = 4000 and L2 = 8000 do 
not seem to differ for (1/6*) = 2000. Fig. |^b shows a case 
of finite size corrections we have dismissed (note the dif- 
ferent scales in the two graphs). It differs from Fig. |l|a 
only by Li = 1000. 

Fig. H illustrates the strong agreement of n{s; 9) at the 
same value of 9 for the same two different sizes L as in 
Fig. |l|a. The two sets of data are virtually indistinguish- 
able, but in this kind of plot it is also almost impossible 
to see a difference between the data of Li — 1000 and 
L2 = 8000, as shown in the inset of Fig. ^ This is also 
the case with the rescaled data below, and the use of 
very large systems throughout this paper might there- 
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fore be "overcautious" in avoiding finite size effects, al- 
though such large sizes are obviously required for an accu- 
rate quantitative analysis of this model. However, when 
it comes only to qualitative analysis, such a judgment 
seems to be justified. On the other hand, an increase 
in system size hardly increases the computing time and 
affects "only" the memory requirements, which forced us 
to implement the algorithm for parallel machines. The 
side effect of using multiple CPUs at the same time is 
a significant reduction of the simulation time especially 
for large values of (1/^), a fact which compensates the 
complications of parallel coding. 

Another indicator for the absence of finite size effects 
is the scaling of the standard deviation of p: If the lattice 
can be split into independent parts, i.e. if subsets of the 
lattice can be considered as independent, the standard 
deviation of p should scale like 1 /L for different values of 
L at given {1/0). Such a behavior can be seen in Tab. |, 
although the standard deviation of p could be calculated 
only roughly. This might explain the slight mismatch for 
{1/0) = 32000, L = 16000,32000. 

For the highest values of (1/6') we could not yet do 
the comparison to another system, so the curve for the 
largest value of {1/9) in Fig. ||is dotted, as their quality 
is not known. However, it is reasonable to assume that 
it is not affected by finite size scaling. 



3.2. The scaling function 

Comparing the different histograms n(s; 9) for differ- 
ent values of {1/6) in a plot enables us not only to find 
the exponent r, but also to find the universal function 
Q as defined in Eqn. |l|. A rough, naive estimate of r is 
given by n{s;9) fitted against s""^, which gives a value 
of r* « 2.1 in our case. Plotting now n{s;d)s'^ double 
logarithmically should allow us to find the "true" value 
of T by performing a data collapse, i.e. choosing r* in 
such a way that horizontal shifts (corresponding to the 
choice of the scale so{9) in the scaling function) make all 
curves collapse. This is shown in Fig. ^, where r* =2.1 
was chosen so that the maxima for the second bumps are 
almost equally high: denoting their position on the ab- 
scissa for each value of 9 by s„ax(^), we have chosen r* 
such that 



n{s^^^{9);9) sC^{9) « const. 



(5) 



According to (|l|) the constant is simply the maximum 
value of Q, namely G{s„^^„{9)/sq{6)), where the value of 
the argument is therefore the same for all 6. 

The value of r* is close to (but not within the error 
of) the exponent found in the literature, r — 2.14(3) 
|, ll (r = 2.15(2) in |l|, r = 2.159(6) in §), which is 
shown in the same figure for comparison. However, it is 
impossible to force the minima (see the down pointing 
marks in Fig. ^ to the same height while maintaining 
the constraint that the maxima remain aligned, i.e. these 
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FIG. 3: The rescaled and binned histogram n(s; S)s^ , where 
T* = 2.10 for (l/e) = 125, 250, 500, • • • , 32000, 64000 (as indi- 
cated) in a double logarithmic plot. The linear size L is cho- 
sen according to the bold printed entries in Tab. | and large 
enough to ensure absence of finite size effects. The error-bars 
are estimated from shorter runs. The rightmost histogram (in 
all figures dotted, (1/9) = 64000) could not be crosschecked 
by another run, see text. Maxima are marked by arrows point- 
ing upwards, minima are marked by arrows pointing down- 
wards. The dashed lines belong to different exponents, whose 
value is specified as the sum of the slope in the diagram and 
T* , i.e. a horizontal line would correspond to an exponent 2.1. 
The shortly dashed line are estimated exponents for different 
regions of the histogram (2.22 within approx. [20, 200] and 
2.19 within [200,2000]), the other exponents are from litera- 
ture, namely 2.14(3) in [|, |] and 223/91 ^ 2.45 in Since 
it was impossible to relate these exponents to any property of 
the data, the exact position of the lines associated with them 
was chosen arbitrarily. 



minima cannot be a feature of the same universal scaling 
function. Otherwise (0) would hold and the quantity 



n{s„,US);9) slJ9) 



(6) 



where s„i„(0) denote the position of the minima, would 
assume the same value for all 9, because they are local 
minima of G, which is supposed to be the same for all 9. 

Since these minima cannot be included in the simple 
scaling defined in (^, they must be explicitly excluded 
by introducing a lower cutoff, so that simple scaling sup- 
posedly sets in only above these cutoffs, excluding espe- 
cially the minima. However, such a lower cutoff would 
apparently have to diverge for {1/9) — > oo - something 
that is certainly beyond any established concept of scal- 
ing. Even when accepting this peculiar scaling behavior, 
a data collapse for the second bump still seems to be 
unsatisfactory, as shown in Fig. ^. 

If one accepts a divergent lower cutoff of the scaling 
function, one has to face the fact that this would describe 
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FIG. 4: The rescaled and binned histogram n{s;9)s'^ , 
versus s/s„ax(^), where t* — 2.10 for (1/S) ~ 
125,250, 500, ••• ,32000,64000 in a double logarithmic plot. 
The scales Srai,^{0) by which the histograms have been shifted 
are the maxima marked in Fig. ^, so that a data collapse 
would be possible. The arrow indicates the order of the data 
in increasing (1/0). 



the behavior of ri in a region, which becomes physically 
less and less interesting in the limit (1/9) oo, because 
the vast majority of events are situated at small s and 
as the second bump moves out to infinity, the scaling 
function hence covers a smaller and smaller part of n. 
However, only a region of ri which covers a non- vanishing 
fraction of events can be physically relevant. 

Concentrating now on the behavior of n up to the min- 
imum (see arrows pointing downwards in Fig. ^), one 
finds that this region is also badly described by a func- 
tion like ([|). First of all, the question of which region 
is supposedly described by the function needs to be an- 
swered. A unique lower cutoff and a dependent upper 
cutoff needs to be found. At first view it looks appealing 
to choose these two marks such that they cover the set 
of data, where the curves fall on top of each other. In 
this case the lower cutoff would be 1 and the upper cut- 
off, 

■Snaivei would have a value smaller than the minima 
marked by downwards pointing arrows in Fig. |^. How- 
ever, this would be described by a function like 

nis; 0) = f{s)g{s/ Snaive) (7) 

rather than (|l|). Note the parameter independent function 
/(s) describing the shape of the curve, while 0(s/sjiaive) 
is a sharp cutoff function. Eqn. does not allow for an 
exponent, /(s) is an arbitrary function. Writing it as 

/(s) = s^'^(ao + higher order corrections) (8) 

defines r to be the steepest descent of this part of the 
curve and gives a value between r^tp = 2.22 and r^itp — 
2.19 (see Fig. |). 



This concept appears to be rather naive - on the other 
hand, it is hard to assume that (||) can still hold: it 
would correspond to (|^) with /(s) replaced by s"^, which 
is a straight line in a double logarithmic plot. Therefore 
(|l|) can apply only to a region in Fig. ^ where the data 
that fall on top of each other form a straight line. Those 
features not already collapsing would then collapse when 
properly tilted (choosing the right r) and shifted (choos- 
ing the right sq). Introducing a lower cutoff at s = 10 
and discarding the data for {1/9) < 2000 then leads to a 
data collapse in a narrow range as shown in Fig. |^. It is 
worthwhile mentioning that even for some 10 < s < 200, 
namely for values of s between the squares and the filled 
circles, none of the data collapse. The exponent used in 
this "collapse" is r^tp = 2.19, as mentioned above. 
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FIG. 5: The rescaled and binned histogram n{s;9)s'^'''^^-, 
versus s/s^i„{9), where Tatp. = 2.19 for (1/0) = 
4000, 8000, 16000, 32000, 64000 in a double logarithmic plot. 
The scales Smi„(0) by which the histograms have been shifted 
are slightly different from the minima marked in Fig. ^, to 
make the collapse as good as possible. The squares and the 
filled circles mark s = 10 and s = 200, respectively, for orien- 
tation and relation to other figures. The arrows indicate the 
order of the data in increasing {1/0). 



By considering the function f{s) it becomes apparent 
that n, and therefore the model, cannot be scale free: 
it depends on the fixed, microscopic scale s = 1. This 
entails that it is always possible to tell (1/6*) by looking 
only at the shape of n; a diagram showing only this shape, 
without any scales on the axes, reveals (1/6*), since a scale 
is intrinsically given by the features of f{s). One would 
only need to rescale and tilt it until it fits the plot Fig. ^ 
and one could identify (1/0). Only if f{s) were scale free, 
i.e. a straight line in a double logarithmic plot, would this 
not be possible. 
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FIG. 6; The position of the minimum in the binned 
and rescaled histogram for different values of r* — 
2.04,2.08,2.10,2.12,2.16. The exponents shown in the plot 
are for comparison only. 




their scaling in (1/6') does not depend strongly on this 
choice. In particular a:„,i„ and x„„x are different for all 
choices of r*. A plot of s^i„{9) versus (1/6') for differ- 
ent values of t* is shown in Fig. For small values of 
(1/0) the minimum is not pronounced enough to survive 
for large values of r* , so these curves do not give a data 
point. Using a linear fit of logs„,i„(6') versus log (1/6') of 
the minimum as found in the rescaled (r*) and binned 
histogram, gives an "exponent" between a;„i„ = 0.93 and 
a^min = 0.98. The same procedure applied to the max- 
ima gives an "exponent" in the range x^^x — 1-18 and 
a^max = 1.22, shown in Fig. ^. One may expect that a;„i„ 
tends towards x^^^ for decreasing r*, as s„i„ increases 
and might enter the scaling region of Smax, but neither 
"exponent" exhibits a systematic variation, and the qual- 
ity of the fit certainly suffers from the rough procedure 
that searches for the extrema in the binned histogram. 
This is unfortunately necessary because of statistical fiuc- 
tuations, in conjunction with the absence of error-bars for 
all data points. 

The scale of the clusters, s„i„/„,ax is related to the cor- 
relation length ^ by the fractal dimension /x, i.e. (see 

i) 

s_/n,.x cx e'^-"'"/-"" . (11) 

Since ^ oc (l/Oy, one should expect v = 
a^min/max/Mmin/max- Thc mlulma arc supposed to be dom- 
inated by smaller, fractal events (see |]l2j), so /i„-,i„ = 
1.96(1) [g and therefore z/„,„ G [0.47,0.50]. The max- 
ima are more likely to be dominated by compact fires, so 
I'max G [0.59,0.61]. It is unclear how the two exponents 
I'min/max 81,16 rclatcd cxactly to the exponents of the two 
correlation lengths found by Honecker and Peschel [^ for 
the connected correlation function i/ = 0.576(3) and for 
the tree-tree correlation function ly = 0.541(4). 

4. DISCUSSION 



FIG. 7: The position of the maximum in the binned 
and rescaled histogram for different values of r* — 
2.04,2.08,2.10,2.12,2.16. The exponents shown in the plot 
are for orientation only. 



3.3. Two length scales 

That n contains features to define at least two scales, 
which apparently diverge in {1/9) with different expo- 
nents, becomes clear when analyzing the scaling of the 
minima and maxima as marked in Fig. ^, using the defi- 
nitions 

s„,je) cx {i/er-^- (9) 

s_x(6') oc (1/6')^— . (10) 

Of course, the exact position of the extrema of n(s; d)s'^ 
depends on its tilt, i.e. on the choice of r*. However, 



Prima facie the DS-FFM looks like a percolation pro- 
cess, and one might naively think that it is indeed a 
percolation process which organizes itself to the criti- 
cal density: sites are occupied randomly and indepen- 
dently and (at least in the thermodynamic limit) there 
is only one cluster which is removed with non-vanishing 
probability, namely the largest. In this way the den- 
sity of occupied sites is automatically reduced below the 
percolation threshold whenever the threshold is reached. 
It is puzzling how remarkably close the tree density in 
the DS-FFM is to the density of empty sites in critical 
site percolation on a square lattice (pffm ~ 0.4078 and 
= 0.40725379(13) ||] respectively). However, the 
removal process involved in thc DS-FFM introduces spa- 
tial correlations which are not present in standard per- 
colation. These correlations are expressed, for example, 
in the form of a patchy tree density distribution [ p2| . 

The purpose of this paper is not to add yet another 
model to the enormous zoo of SOC models. However, in 
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FIG. 8: The rescaled and binned histogram n{s;6)s^ 
(again r* — 2.10), for a modified model, where the largest 
cluster in the system is removed after each driving step, 
for (l/e) = 1000,2000,4000 (as indicated) with linear sizes 
L = 2000,2000,4000. The inset shows the same data on 
the scale of Fig. ^ for comparison. The data for (1/9) — 
1000, 2000, 4000 of the original model as shown in Fig. ^ are 
dotted. The pecuhar behavior of the different height- scaling 
of the minimum and the maximum is again visible. 



order to investigate certain features of the given model 
and identify underlying mechanisms, it makes sense to 
modify it slightly. The outcome for the histogram of 
the DS-FFM modified such that the largest cluster is re- 
moved after each driving step is shown for a few values 
of (1/6') in Fig. |[ The distinctive feature of a mini- 
mum which scales differently from the maximum is again 
present, as the peaks of the maxima have approximately 
the same height, while the height of the local minima 
varies among different values of 9. The inset of this figure 
shows the histogram on the same scale as Fig. ^ together 
with the data of the original model (dotted) with the cor- 
responding values of 9. One can understand that they do 
not fall on top of each other, because the relaxation rule 
in the modified model erases much larger clusters than 
in the original model. 

Fig. ^ shows a second modification of the model, where 
again the largest cluster is removed during relaxation and 
in addition the driving is changed such that the density 
of trees, p, is the same before each relaxation; the trees 
removed during the relaxation are just filled in randomly 
afterwards. This model differs from standard percola- 
tion only by its updating scheme [20[ |. In order to com- 
pare the outcome with the original model, the values of 
p have been chosen close to the values given in Tab. |. 
Indeed, the feature of different scaling of the extrema is 
still present, but it disappears completely if the density 
is increased to Ppcrc = 0.592746 which is shown in 
the same figure as the large bump. This curve does not 



FIG. 9: The rescaled and binned histogram n(s; p)s^ (again 
r* — 2.10), for a modified model, where the largest cluster 
in the system is removed in each relaxation step and the 
corresponding number of trees is filled back into the sys- 
tem afterwards. The three small values of p chosen, p = 
0.3975, 0.4005, 0.4025 correspond (up to the last digit) to the 
values of the tree density for (1/61) = 1000, 2000, 4000 respec- 
tively, see Tab. |. The linear size was L = 1000, 2000, 4000. 
The corresponding data of the original model are shown dot- 
ted. The peculiar behavior of the different height-scaling of 
the minimum and the maximum is again visible (a correct 
tilt T* would make it even more pronounced), but disappears 
obviously for p = Pperc - for these data it is relevant to men- 
tion that n{s) was measured after the relaxation. The filled 
circles show the exact results for the lattice animals |^ |l^, p^ ] 

at p = Ppcrc. 



vary much if a much smaller system size is simulated at 
this density, so we expect it essentially to be free of finite 
size corrections. Since it represents a correlated percola- 
tion process, it is just consistent that this bump does not 
cover the exact results for the lattice animals of standard 
percolation [|[ ^ ^ at p = Ppcrc shown as filled circles 
in Fig. ^. The dotted graphs in the figure show the corre- 
sponding data of the original model. Again they do not 
match apart from the region of very small s. Unfortu- 
nately the simulations of the so-modified model are very 
expensive in CPU time, because the mass of the largest 
cluster needs to be refilled after each relaxation, so that 
only 50.000 updates for transient and statistics could be 
done. 

Since the feature of different scaling survives the mod- 
ifications described above, it seems reasonable to assume 
that any relaxation rule that favors the largest cluster 
leads to the peculiar behavior. Its disappearance at high 
densities can be explained by the extremely small cut- 
off in the distribution, which leads to a domination of 
the statistics by very small clusters, while a single, enor- 
mously large one dominates the burning (the average size 
of the burnt cluster for p = pporc was 355811). However, 
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much more careful and detailed investigations of mod- 
els like the modification described above are required to 
gain a full understanding of the underlying mechanisms. 
In particular, this should include a modification of the 
rules such that the feature disappears. 

Honecker and Peschel have calculated the correla- 
tion length not only for the probability that two sites be- 
long to the same cluster, but also for the probability that 
two sites are occupied at all. The correlation function for 
the latter is of course a i5 peak in ordinary percolation, as 
there are no spatial correlations for the distribution of oc- 
cupied sites by construction. However, in the DS-FFM 
the correlation length for this quantity, ^, is finite and 
seems to diverges when approaching the critical point. 
It is highly remarkable that Honecker and Peschel con- 
clude from their simulations that this correlation length 
diverges slightly slower than the correlation length of 
the probability for two sites to belong to the same clus- 
ter. This seems to indicate that for sufhciently large 
scales the spatial correlation of the occupation proba- 
bility becomes arbitrarily small, so that on sufficiently 
large scales the DS-FFM occupation is uncorrelated and 
therefore tends to standard percolation. In other words, 
it seems to be possible to rescale or "renormalize" the 
DS-FFM to make the occupation correlation arbitrarily 
small. This would introduce higher order interactions, as 
known from standard real space renormalization group 
and would explain the difference in critical density be- 
tween the rescaled DS-FFM and standard percolation. 
However, if this "mapping" is valid, one should find the 
exponent for the divergence of ^s/S, to be the same as in 
standard percolation, but this is precluded by numerics. 

It has been suggested at least twice |l2| , that the DS- 
FFM is a superposition of cluster distributions np„„(<, p) 
of standard percolation for a whole range of concentra- 
tions p, weighted by a certain distribution function w{p), 
i.e. dpw{p)n{s,p). Obviously such an assumption ne- 
glects spatial correlations. We recall the following result 
from standard percolation theory 

n{s,p)^s-^C{^s/{p^p,)-'/^) , (12) 

where C denotes the cutoff function and the exponents a 
and Tperc have their standard definitions. Assuming that 
the weighting function w{p) is analytic around the critical 
concentration in standard percolation, pc, ( |l^ ) leads to 

[ dpw{p)n{s,p) (X s-^''^''"'+'''> . (13) 
Jo 



This gives rise to an exponent r = 223/91 ~ 2.45, how- 
ever, this is definitely not supported by numerics (see 
Fig. I). 

It remains completely unclear how to characterize the 
scaling of the DS-FFM in two dimensions. Apparently it 
is not a mere superposition of two simple scalings, as re- 
cently speculated ||l2| . Moreover the model does not seem 
to be scale free as described above and it does not seem 
to be possible to identify a unique power law behavior of 
the cluster size distribution. Nevertheless effective power 
law behavior over restricted regions has clearly been pro- 
duced by the model, making it potentially relevant to 
observation. 

All we can conclude is that the DS-FFM is not critical 
in the sense of simple scaling. It reminds us that a diver- 
gent moment (here (s), the second moment) can be re- 
garded as a unique sign of emergent scale invariance only 
if we are certain that one single scale is sufficient to char- 
acterize the system. If there is more than one relevant 
scale, different properties of the system might depend on 
different scales which may or may not diverge. 
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